function [corrCell]=correlation_AVA(type, seisReal, seisSim, layers)
%% COMPUTE THE CORR COEF BETWEEN ALL THE SIMULATIONS WITH THE REAL SEISMIC
%
% type - TYPE OF CORRELATION COMPUTED

% seis_real [cell] - REAL SEISMIC
% seis_sim [cell]  - SIMULATED SEISMIC (just the data from the "noIt" iteration)
% layers           - NUMBER AND SIZE O EACH LAYER USED FOR THE ITERATION
%
% corr_cell [cell] - [size(layes,1) noAngles]
%
% [MODIFICATIONS]
%   AUG 2019    - Major performance update. Remove a loop on layers.
%   JUL 2019    - CAP ON cc IS DONE OUTSIDE. Handle NullDataValues
%   04 DEC 2015 - MAJOR UPDATE ON PERFORMANCE - got ride of for loops
%   21 APR 2014 - CORRELATION VOLUME AS ONLY size(layers,1) SIZE
%   29 SEP 2013 - BUG FIX. CONSIDER ONLY THE SYNTHETIC AT CURRENT SIMULATION
%   10 JUL 2013 - BUG FIX. WHEN REAL SEISMIC == 0 WITHIN LAYER ind WAS NOT
%                   UPDATING. PERFORMANCE UPDATE
%   05 MAR 2013 - NEW CONDITION TO CHECK IF ORIGINAL SEISMIC IS == 0 IN THE
%                   ENTIRE LAYER
%   03 OCT 2012 - LAYERS ARE NOW COMPUTED OUTSIDE IN MAIN
%   11 SEP 2012 - GLOBAL CORRELATION IS COMPUTED CORRECTLY IN A NEW FUNCTION
%   07 SEP 2012 - PEARSON CORRELATION COMPUTED BY HAND
%   04 SEP 2012 - CHANGES FOR SINGLE SIMULATIONS RUN - MEMORY UPDATE.
%   21 AUG 2012 - MEMORY ALLOCATION UPDATE
%   02 AUG 2012 - COMPUTE AND WRITE INTO LOG FILE THE GLOBAL CORRELATION OF
%                   SYNTHETIC PER SIMULATION
%   20 JUL 2012 - DEFAULT PARAMETERS FOR MIN AND MAX LAYER SIZE REMOVED;
%                   CORRELATION VALUES ARE ALWAYS BTW. 0 AND 1
% 11 JUL 2012 - TYPE ADDED TO RUN WITH QUASI CORRELATION
% 21 JUN 2012 - QUASI CORRELATION ADDED
% 27 APR 2012 - FIXED LAST LAYER IS ALWAYS > MinLayerSize
% 23 APR 2012 - NOW COMPUTES CORR COEF BTW LAYERS WITH VARIABLE SIZE
% 30 MAR 2012 - NOW COMPUTES CORR COEF BTW LAYERS WITH FIXED SIZE
%
% LA - 22 APR 2012
%

%% allocate and initialize
sizeGrid    = size(seisReal);
sizeGather  = size(seisReal{1,1});

seisReal    = cell2mat(seisReal);
seisSim     = cell2mat(seisSim);
corrMat     = zeros(size(layers,1) * sizeGrid(1), size(seisReal,2));
layersAll   = repmat(layers, sizeGrid(1), 1);

%% COMPUTE CORRELATION BTW LAYERS

if type == 0 % CORRELATION COEFFICIENT
    % index for the layering
    indLayer    = 1;
    for n=1:size(layersAll,1)
        seisRealLayer   = seisReal(indLayer: indLayer+layersAll(n)-1,:);
        seisSyLayer     = seisSim(indLayer:indLayer+layersAll(n)-1,:); 
        if (sum(seisRealLayer) == 0 & sum(seisSyLayer) == 0)
            corrMat(n ,:) = 1;
        else
            corrMat(n ,:) = (sum((seisRealLayer - repmat(mean(seisRealLayer),size(seisRealLayer,1),1))...
                .*(seisSyLayer - repmat(mean(seisSyLayer),size(seisSyLayer,1),1))))./...
                (sqrt(sum((seisRealLayer - repmat(mean(seisRealLayer),size(seisRealLayer,1),1)).^2))...
                .*sqrt((sum((seisSyLayer - repmat(mean(seisSyLayer),size(seisSyLayer,1),1)).^2))));
        end
        indLayer=indLayer+layersAll(n,1);
    end
    
elseif type == 1

    indLayer    = 1;

    for n=1:size(layersAll,1)
        seisRealLayer   = seisReal(indLayer: indLayer+layersAll(n)-1,:);
        seisSyLayer     = seisSim(indLayer:indLayer+layersAll(n)-1,:); 

        corrMat(n,:) = (2*sum((seisRealLayer) .*seisSyLayer,'omitnan')) ...
            ./ (sum(seisRealLayer.^2,'omitnan')...
            + sum(seisSyLayer.^2,'omitnan'));

        indLayer=indLayer+layersAll(n,1);
    end
end

corrMat(corrMat < 0) = 0;
corrMat(isnan(corrMat(:,:)))=0;
corrCell = mat2cell(corrMat, repmat(size(layers,1),sizeGrid(1),1),  repmat(sizeGather(2),sizeGrid(2),1));

end
